Ab initio investigation of competing instabilities in ferroelectric perovskite PbTiOa 
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We have developed first-principles models, based on a general parametrization of the full potential- 
energy surface, to investigate the lattice-dynamical properties of perovskite oxides. We discuss 
the application of our method to prototypic ferroelectric PbTi03, showing that its properties are 
drastically affected by a competition between structural instabilities. Indeed, we confirm previous 
indications that the destructive interaction between the ferroelectric and antiferrodistortive (involv- 
ing rotations of the C>6 octahedra) soft modes shifts PbTiCVs Curie temperature by hundreds of 
degrees. Our theory provides unique insight into this gigantic effect and its dynamical character. 



PACS numbers: 63.70.+h,64.60.De,77.80.B-,71.15.Mb 

Perovskite oxides are a large family of materials of 
great fundamental and applied interest. In many cases, 
the structural and lattice-dynamical features of the com- 
pounds are critical to determine their properties. Most 
notably, this is the case of ferroelectric (FE) perovskites 
- whose spontaneous polarization is usually the result of 
a structural phase transition [l a - and related compounds 
such as magnetoelectric multiferroics - whose properties 
can be greatly enhanced by engineering the lattice re- 
sponse to external fields [5J. First-principles studies of 
lattice-dynamical phenomena are usually restricted to 
the limit of low temperatures, as the inclusion of ther- 
mal effects requires large simulation boxes (thousands of 
atoms to get a realistic description of most properties of 
interest) and a good sampling of configuration space (tens 
of thousands of Monte Carlo sweeps to compute statisti- 
cal averages). Hence, it is not yet feasible to tackle such 
problems from first-principles directly. 

The situation greatly improved with the development 
of effective models that (1) capture in a general and 
mathematically simple form the energy surface associ- 
ated with the most relevant (i.e., lowest-energy) struc- 
tural distortions and (2) whose parameters are com- 
puted from density-functional theory (DFT). This so- 
called effective Hamiltonian approach has been success- 
fully applied to FE crystals like BaTi0 3 [3] and PbTi0 3 
(PTO) @], solid solutions like PbZr!_ x Ti x 3 (PZT) 
[51 [3], and even multiferroics like BiFe0 3 [7J. However, 
by neglecting all but the lowest-energy degrees of free- 
dom, this scheme may introduce quantitative (e.g., in 
the description of thermal expansion [8j) and qualitative 
(e.g., whenever secondary modes are important [9]) er- 
rors. Further, the approach does not seem well suited 
to tackle structurally complex cases - from domain walls 
to heterostructures displaying interface-driven phenom- 
ena |10| - in which it may be difficult to identify a 
small set of distortions that capture the main effects. 
Shell models [TT] and other schemes [12] providing a full 
atomistic description have also been derived from first- 



principles for these materials; however, such approaches 
rely on specific forms of the interatomic potentials, which 
may restrict their accuracy and generality. 

To overcome these limitations, we have developed an 
aproach that provides a full atomistic description of the 
materials while retaining the general energy parametriza- 
tion of the effective-Hamiltonian scheme. Here we de- 
scribe our first application - to the study of PTO's FE 
phase transition -, which illustrates the predictive power 
and unique insights that the new method offers. 

Methods.- Let us briefly describe our all-atom effec- 
tive Hamiltonians, which will be presented in detail else- 
where [13 j . As in the usual approach [3], we write the 
energy of the material as a Taylor series around a (cu- 
bic perovskite) reference structure. Our Taylor series, 
though, is a function of all the atomic degrees of free- 
dom, not only the lowest-energy ones. It is convenient to 
split our variables into atomic displacements (generically 
denoted by u) and global strains of the simulation box 
(77), so that we can write the energy as 

E(u, v)=E + E(u 2 ) + E{jf) + E(r,u 2 ) 

+ E(u 3 ,u\...) + E(r] 3 ,r, 4 ,...) (1) 
+ E(r l 2 u 2 ,r l 3 u 2 ,...,r i u 3 ,...). 

The first line in Eq. (1) contains the harmonic terms in 
u (i.e., the force- const ant matrix K) and 77 (frozen- ion 
elastic tensor), which can be readily obtained from most 
first-principles codes; we computed them for PTO using 
density-functional perturbation theory as implemented in 
ABINIT [MlflT] . Figure 1 shows the dispersion curves for 
the K qs eigenvalues of K, which include a variety of struc- 
tural instabilities of the cubic phase (i.e., distortions with 
Kqs < 0) and are exactly captured by our model [15]. The 
first line of Eq. (1) also includes the lowest-order non- 
zero rj-u couplings [191 . which we obtain from the force- 
constant matrices computed for slightly strained cells. 
Such a strain-phonon coupling term is the only rj-u inter- 
action considered in the usual effective Hamiltonians of 
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FIG. 1: (Color online.) Dispersion curves for the K qa eigen- 
values of the force-constant matrix K of PTO's cubic phase. 
Lines and circles are the results obtained ab initio and with 
our effective model, respectively. The line color indicates the 
atomic character of the eigenmode; red, green, and blue cor- 
respond to Pb, Ti, and O, respectively. Unstable (n qs < 0) 
FE and AFD modes are sketched. The AFD a (AFDi) modes 
involve an antiphase (inphase) modulation of the Og-rotations 
along the rotation axis. 

ferroelectrics. 

The second and third lines of Eq. (1) include anhar- 
monic interactions that in principle can extend to ar- 
bitrary order and need to be truncated [IS]. We com- 
puted the parameters in such terms by fitting to the 
key features (energy, structure, X-eigenmodes at se- 
lected q points) of some relevant low-symmetry low- 
energy phases, which include several structures in which 
FE and antiferrodistortive (AFD) modes are condensed. 
[As described in Fig. 1, AFD modes may involve inphase 
(AFDj) and antiphase (AFD a ) rotations of the 0§ octa- 
hedra; only the AFD a modes are relevant in the present 
discussion.] We first considered low-symmetry phases ob- 
tained by imposing the cubic cell (rj — 0) to compute 
the parameters in the E(u 3 , it 4 , ...) term; we concluded 
it can be well approximated by (1) truncating the series 
at fourth order and (2) including only pairwise interac- 
tions between nearest-neighboring atoms (Pb-0 and Ti- 
O). We then considered the fully-relaxed low-symmetry 
structures, and found that high-order rj-u couplings are 
needed in order to describe their main features accurately 
and, simultaneously, reproduce the first-order character 
of PTO's FE transition. (Our model includes selected r\-u 
couplings up to 7/ 4 u 2 order.) The E(r) , rj A , ...) term, on 
the other hand, could be neglected. The highest-order 
terms were chosen to be even powers of u and ry, and 
the corresponding parameters restricted to be positive, 
to guarantee the global stability of the model. 

Typically, the effective model of PTO used in this 
work reproduces the first-principles energies of the rel- 



evant low-symmetry phases within a 15%, and their low- 
lying Jf-eigenvalues within 0.5 eV/A 2 . Note that our 
parameter-fitting procedure implicitly captures the usual 
features of atomistic models of ferroelectrics. For ex- 
ample, the double-well potential associated with the FE 
instability stems from (1) the unstable polar modes in- 
cluded at the harmonic level and (2) the anharmonic 
terms fitted to describe the stable low-symmetry phases 
obtained ab initio. The anharmonic interactions between 
modes are captured by fitting to the n qs spectrum of dis- 
torted structures, etc. 

Basic results.- We solved our model by performing 
Monte Carlo simulations, typically using 20000 sweeps 
for thermalization and 80000 for computing averages 
at each considered temperature (T). The simula- 
tion box contained 8x8x8 perovskite cells (i.e., 2560 
atoms) and periodic boundary conditions were employed. 
(10x10x10 supercells of 5000 atoms were used for T ss 
Tc, where size effects become more relevant.) Figure [5] 
shows our basic results. Our model predicts that, at a 
temperature Tq ~ 325 K, PTO undergoes a sharp transi- 
tion between the high-T cubic (Pm3m) and low-T tetrag- 
onal (Pimm) phases. This transition carries with it the 
deformation of the unit cell [Fig. 2(a)] and the develop- 
ment of an spontaneous polarization P [Fig. 2(b)]. 

We confirmed that the polarization is the primary or- 
der parameter of the transformation by running simula- 
tions with the cubic cell fixed (77 = 0); in that case we still 
obtain a FE transition at a similar Tq [inset of Fig. 2(b)], 
but the polar phase has now a rhombohedral symmetry 
(i?3m) with P x = P y = P z . Further, at T w 200 K 
we observe a second transition in which rotations of the 
6 octahedra around the polar axis are condensed (the 
symmetry is i?3c and the cell doubles). These results 
show the importance that the rj-u couplings have in de- 
termining PTO's behavior, and suggest that mechanical 
constraints (e.g., epitaxial strain) can greatly affect it. 

Figure 2(c) shows the T-dependence of the atomic po- 
sitions in PTO's unit cell. We find that, while the magni- 
tude of the FE distortion depends strongly on T, the dis- 
placement pattern remains relatively constant. Hence, in 
this sense our results validate the usual assumption that 
the main features of the FE transition can be captured by 
effective models including only one type of polar distor- 
tion. Note also that at T rj Tc the atomic displacements 
become soft [i.e., there is a large spread of their instan- 
taneous values, as shown in the inset of Fig. 2(c)], which 
leads to an enhancement of the dielectric and (below Tc) 
piezoelectric responses (not shown here). 

Because it includes all the atoms, our model should 
correct the usual effective-Hamiltonian underestimation 
of the thermal expansion [8 . Indeed, above Tq we obtain 
an expansion coefficient a = 18.2 x 10~ 6 °C _1 , which 
exceeds the experimental value of 12.6 x 10 -6 °C _1 [20j . 

What controls Tq ?- Our computed Tc is about 325 K, 
while the experimental result is around 760 K |20j . This 
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FIG. 2: (Color online.) T-dependence of the quantities characterizing PTO's FE transition. Panel (a): strains and volume 



change with respect to the reference cubic cell. Panel (b): Spontaneous polarization P a = v~ L Ylip z %p : aUip, where v is the unit 
cell volume, a and /3 are spatial directions, i labels atoms within the unit cell, Zip, a are Born charges, and Uip are supercell- 
averaged atomic displacements. Inset: results of simulations at fixed n — 0; here, the FE and AFD a distortions are obtained by 
projecting the thermal-averaged configuration onto the corresponding unstable if-eigenmodes of the cubic phase (see Fig. [TJ ; 
units are arbitrary. Panel (c): Supercell-averaged displacements Ui a along the direction of the spontaneous polarization. A 
global displacement is chosen so that J^. u ia = 0. Inset: standard deviation of the Pb and Ti displacements, in arbitrary units. 



severe underestimation is surprising, as previous stud- 
ies of PTO based on less-complete theories had shown a 
much better agreement with experiment. In particular, 
Waghmare and Rabe (WR) constructed a model that 
neglects all degrees of freedom but the soft polar modes 
(treated as lattice Wannier functions) and cell strains, 
and obtained a value of 660 K [21] . 

We checked our model reproduces the energetics of 
the FE instabilities given by the WR Hamiltonian quite 
closely, despite the differences (e.g., DFT functional, 
pseudopotentials) in the ab initio calculations employed 
to compute the parameters. Further, we ran simulations 
with modified versions of our model to test subtle features 
of the WR energy parametrization (e.g., the inclusion of 
high-order terms for the FE distortions), and concluded 
that they have a negligible effect on Tq. 

We thus turned our attention to the qualitatively dis- 
tinct features of our model. Most notably, we describe 
not only the FE instabilities, but also the unstable AFD 
distortions shown in Fig. 1. It is known that, in most 
perovskite oxides, the interaction between FE and AFD 
modes is a competitive one, so that they tend to sup- 
press each other. The best studied case may be that 
of SrTi0 3 (STO), for which Zhong and Vanderbilt pre- 
dicted, by means of an effective Hamiltonian approach, 
that the temperature of STO's AFD-related structural 
transition would be about 25% higher if the soft FE dis- 
tortions of the material were suppressed [22] . 

The histograms in Fig. [3] show that similar effects 
occur in our PTO simulations. Above Tq, the AFD 
modes have relatively large amplitudes, comparable in 



magnitude with those of the FE distortions. In con- 
trast, they get significantly suppressed below Tq, in a way 
that clearly reflects the breaking of the cubic symmetry: 
the z-oriented spontaneous polarization restrains more 
strongly the transversal 06-rotational modes AFD ax and 
AFD aj( . Note that the simultaneous occurrence of FE a 
and AFD a ^ distortions, with a ^= /3, results in a signifi- 
cant deformation of the Oq octahedra, which may explain 
why these transversal FE-AFD interactions are particu- 
larly destructive. 

To investigate the effect of the AFD distortions on the 
FE transition temperature, we ran simulations in which 
the 06 rotations were cither totally (no- AFD case) or 
partly (only-AFD z case) suppressed. We imposed these 
constraints by restricting the motion of the oxygen atoms 
as shown in the sketches in Fig. [4] Let us stress that 
these constraints do not affect the energetics associated 
with the development of the spontaneous polarization, 
the FE ground state being exactly retained. (Even in 
the only-AFD z case - in which 06 rotations are allowed 
only around the z axis, which breaks the cubic symme- 
try of the model -, we still keep the six equivalent FE 
minima with P along the ±a;, ±y, and ±z directions.) 
Figure [4] shows the results: In the no- AFD case we ob- 
tain Tq w 650 K, which is very close to the WR result. 
We also ran this case at fixed n = (inset of Fig. |4| and 
obtained a similarly high Tq. In the only-AFD z case the 
calculated Tq is about 500 K, and the spontaneous polar- 
ization is found to point specifically along z. (We checked 
this by running several independent simulations starting 
from the u = r\ = configuration.) Note that these effects 
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FIG. 3: (Color online.) Histograms of the FE and AFD a dis- 
tortions, which are quantified by projecting the instantaneous 
atomic structures on the corresponding unstable eigenmodes 
of the cubic reference phase. 

- i.e., the shift in Tq and the preferential P-direction in 
the only-AFD 2 case - have a strictly dynamical origin. 
For example, in our only-AFD z simulations, the AFD a2 
component fluctuates around its zero average value and 
hampers the development of the polarization, especially 
along the perpendicular directions x and y. Thus, P z is 
dynamically favored and condenses at To 

A couple of additional points are worth making. (1) 
We find that strain does not play a big role in how the 
FE-AFD competition affects Tq , as our regular and fixed- 
cell simulations lead to very similar results in this regard. 
Note that, in principle, one could have expected other- 
wise: Because the FE and AFD instabilities tend to re- 
act differently to an external pressure [23l EH |26] , a FE 
distortion should imply cell deformations that, in turn, 
should suppress the AFD modes, and viceversa. How- 
ever, our simulations suggest that such an effect is not 
very important, probably because at the temperatures at 
which the FE-AFD competition is relevant - i.e., in the 
range from 350 K to 650 K where the material would be 
FE in absence of AFD modes -, the correlations between 
atomic displacements and strains are relatively small. 

(2) Our claim that a destructive FE-AFD interaction 
exists in the fixed-cell case may seem surprising, as a com- 
bined FE+AFD ground state occurs for -q = [see inset 
of Fig. 2(b)]. However, note that competing instabilities 
may not fully suppress each other, and that is the case 
here. To make the argument more quantitative, let us 
note that the fixed-cell FE phase lies —59 meV/f.u. be- 
low the cubic reference structure, while we get the AFD 
phase at —56 meV/f.u. The combined FE+AFD struc- 
ture has an energy of —62 meV/f.u., far above the value 
of —115 meV/f.u. (where 115 = 59+56) that one would 
expect if the FE and AFD modes did not interact at 
all [23]. It is thus clear that these two instabilities com- 
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FIG. 4: (Color online.) T-dependence of the polarization 
(only one component is shown) as obtained in the regu- 
lar (solid squares), only-AFD 2 (open squares), and no- AFD 
(open circles) cases. The arrows in the sketches indicate the 
oxygen displacements that are allowed in the only-AFD z and 
no- AFD cases; the motion of the Pb and Ti atoms is uncon- 
strained. Inset: Results from simulations with fixed n = 0; we 
show the spontaneous polarization distortion for the regular 
(filled squares) and no-AFD (open circles) cases. 



pete in the fixed-cell case, even if they coexist at low 
temperatures. [In essence, in the regular case the polar 
distortion is larger (1.01 C/m 2 at K, as compared with 
0.57 C/m 2 when rj = 0) and the FE-AFD interaction fully 
suppresses the AFD instability] 

Let us note that a very large competition-driven shift 
in Tq had already been predicted (but not much empha- 
sized) for PZT by Kornev et al. [Jjj. The effective model 
used by these authors included some particular FE-AFD 
couplings that allowed them to rectify the very large Tq 
obtained previously from a simpler theory [5]. In con- 
trast, our Hamiltonian captures the FE-AFD couplings 
in an implicit and non-specific way; in fact, when con- 
structing our model we were not aware of the importance 
of this interaction, nor did we need to assume any form 
for it. Hence, our work ratifies that gigantic competition- 
driven effects can occur in PTO-related compounds. In- 
terestingly, the analogous phenomena in SrTi03 , which is 
generally considered the prototypic example for compet- 
ing structural instabilities, are probably much smaller; in 
particular, Zhong and Vanderbilt obtained shifts in the 
transition temperatures of about 35 K |22j . 

The origin of the large difference between our com- 
puted Tc (325 K) and the experimental one (760 K) 
remains an open question. We ran additional tests to 
check the effect of some of the approximations made in 
our model [13] . and did not observe any significant im- 
provements. Hence, we tend to attribute the disagree- 
ment to the first-principles methods used to compute our 
Hamiltonian parameters. Indeed, to obtain the right Tq, 
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the ab initio simulations should describe accurately (1) 
the strength of the FE instabilities, (2) the strength of 
the AFD instabilities, and (3) the magnitude of the FE- 
AFD coupling. Our results suggest that quantity (1) may 
be underestimated, and (2) and (3) overestimated, by 
our current calculations, despite our using the best DFT 
methods available for this task (i.e., recent functionals 
that have been shown to perform very well for these ma- 
terials [23 [2H])- Clearly, T c 's large sensitivity to the 
FE-AFD competition sets a very stringent requirement 
for the accuracy of the ab initio calculations. Thus, our 
work suggests that, in spite of recent progress, we still 
lack DFT methods that describe accurately the thermo- 
dynamic properties of FE materials like PTO. 

Summary.- We have introduced a novel approach 
for first-principles investigations of the lattice-dynamical 
properties of perovskite oxides. More precisely, we have 
developed effective models that are based on a general 
parametrization of the energy of the material (i.e., we 
make no a priori assumptions on the form of the inter- 
atomic couplings) and include all atomic degrees of free- 
dom. The application of the new scheme to PbTiC>3 has 
allowed us to investigate the competition of instabilities 
at play in this compound, which turns out to have dra- 
matic effects in its properties. Indeed, we found that 
such a competition results in a reduction of PbTiCVs 
Curie temperature by as much as 300 K. Our simula- 
tions provide unique insights into this gigantic effect and 
its dynamical character. 
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